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Abstract 



We investigate diffusion-limited reactions between a diffusing particle and a target site on 
a semiffexible polymer, a key component of DNA-protein interaction and polymerization 
of cytoskeletal filaments. Our theory focuses on two competing effects: internal polymer 
fluctuations, which speed up association, and hydrodynamic coupling between the particle 
and chain, which inhibits association. Internal polymer motion is described using a mean 
field dynamical theory, while the hydrodynamic coupling between polymer and particle is 
incorporated through a simple heuristic approximation. Both of these we validate through 
comparison with Brownian dynamics simulations. Neither of the effects has been fully con- 
sidered before in the biophysical context, and we show they are necessary to form accurate 
estimates of reaction processes. The association rate depends on the stiffness of the poly- 
mer and the particle size, exhibiting a maximum as a function of persistence length and a 
minimum with particle radius. In the parameter range relevant to DNA-protein binding, the 
rate increase is up to 100% compared to the Smoluchowski result for simple center-of-mass 
motion. The quantitative predictions made by the theory can be tested in single molecule 
experiments. 
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1 Introduction 

Reactions between semiflexible polymers and small molecules are ubiquitous in cells, playing 
a crucial role in a large number of biological processes: examples include the interaction of 
gene-regulating proteins with specific target sites on DNA, and the polymerization of DNA or 
structural proteins such as actin and tubulin. Many of these reactions are diffusion-limited — 
the activation energy is negligible compared to the thermal energy k^T — so the overall 
association speed therefore depends on the rate at which the reactive molecules approach 
each other. 

DNA-protein interaction has been the most widely studied process of this type (0,0), at- 
tracting attention since the first measurement of the reaction rate between the Lac repressor 
and operator (J3|) revealed that it far exceeds the 3D diffusion limit. The quest to identify 
the underlying mechanism culminated in the seminal idea of facilitated diffusion by Berg, 
Winter, and von Hippel (jl, 0, @), but the general description of DNA-protein interaction 
is still far from complete. Although recent single molecule experiments show evidence that 
certain proteins indeed make use of facilitated diffusion (0, @), the majority of measured 
reaction rates for DNA-binding proteins are on the order of the Smoluchowski result for 3D 
diffusion (Q], [qI). Thus recent years have seen extensive theoretical efforts (!(], U, 12, 13, 14) 



revisiting the underlying assumptions of facilitated diffusion, and examining it anew in re- 
sponse to experimental advances. 

In the current work we focus on two aspects that enter into any theory of facilitated 
diffusion (or polymerization): for a bimolecular reaction involving a specific site on the 
semiflexible polymer and a free particle, how does the association rate depend on (i) the 
conformational fluctuations of the polymer contour in solution; (ii) the long-range hydrody- 
namic coupling between the polymer and the particle. Understanding these two intimately 
linked effects is necessary for any quantitative prediction of empirical rates, but they have 
been largely neglected in the context of DNA-protein association. Diffusion-controlled re- 
actions with targets on flexible polymers are a well-established subject in polymer physics, 



with general analytical frameworks developed by Wilemski and Fixman (Il5l . Il6l ). and Sz- 



abo, Schulten, and Schulten (1171 ) decades ago. Using the Wilemski-Fixman approach and a 



Gaussian model for circular polymers, a study by Berg looked at the influence of internal 



DNA motion on the binding of proteins (18). But little is known about the role of chain 



stiffness, and both Berg's study and the classic polymer reaction rate theories do not in- 
clude hydrodynamics. This is a significant oversight, because solvent-mediated interactions 
modify not only the approach of the particle to the target site, but also the fluctuation of 
the entire polymer contour, and thus all time scales involved in the dynamics. Moreover, 
the two effects are in opposition: polymer fluctuations lead to the target site exploring a 
larger configurational space than simple center-of-mass diffusion, and hence the association 
rate is enhanced; the hydrodynamic coupling between the particle and coil, on the other 
hand, reduces their relative mobility, decreasing the rate. The observed binding rate is a 
subtle competition between these two phenomena. As experiments are providing an ever 
more detailed picture of biological reactions at the single molecule level, we need theories 
that begin to grapple with the full complexity of polymer-particle diffusive motion. 

The paper is organized as follows: in Sec. [2] we briefly describe the Brownian dynamics 
(BD) simulation method; this is used to independently test the theoretical description of 
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polymer-particle dynamics in Sec. [3] The latter is divided into three parts: Sec. [3J] presents 
a short review of the mean field theory for semiflexible polymers, which gives an accurate 
description of the internal fluctuations of the chain; Sec. 13.21 considers the hydrodynamic 
coupling between the free particle and the polymer coil; finally Sec. 13.31 describes the renewal 
approach used to derive first passage times and association rates. The detailed validation 
of the theory through BD simulation results is presented in Sec. HI We then investigate 
how the mechanical characteristics of the semiflexible polymer — the contour and persistence 
lengths — affect the reaction rate. Ultimately our theory yields quantitative predictions that 
can be tested empirically, and Sec. [5] summarizes the main results and offers suggestions for 
future experiments. 



2 Brownian dynamics simulations 



To simulate a polymer in solution we adopt a standard Brownian dynamics (BD) scheme (1191 ). 
modeling the polymer as a chain of M beads of radius a. For the low Reynolds number 
regime, the Langevin equation governing the time evolution of the position V{ of bead i is 
given by 



drjjt) 
dt 



M 
3=1 



dU(n 



The long-range hydrodynamic interactions — the fact that a force fj acting on bead j creates 
a flow- field affecting the motion of bead i — are described by the Rotne-Prager mobility matrix 



fMoSij 1 + (1 - 5ij)*jt(rij), 



1 + 



+ 



1 

y 



y 



*3 
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where = — rj, = \rij\, 1 is the 3x3 identity matrix, and /io = (67r?7a) _1 is 
the Stokes self-mobility of a sphere of radius a in a solvent of viscosity rj. The stochastic 
contributions £j(t) to the equation of motion are assumed to be Gaussian random vectors 
which are hydrodynamically correlated according to the fluctuation-dissipation theorem: 



(4) 



The inter-bead potential U = Uwlc + Ulj determining the configuration-dependent forces, 
where 
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consists of a shifted harmonic potential between adjacent beads of strength 7 = 200k-sT/a, a 
bending potential of strength k = l v k^T between adjacent bonds, and a pairwise truncated 
Lennard- Jones potential C/lj of strength w = 3k^T. Here 6i is the angle between the bond 
vectors and r i+ i ti . The first term in the worm-like chain potential tAvLC keeps the 

contour length L approximately fixed, while the second one with modulus k takes care of the 
bending stiffness of the chain. The repulsive Lennard- Jones potential L/lj prevents significant 
bead overlap, which is a source of numerical instabilities. In Sec. 14.31 where we model the 
motion of a free particle (i.e. a protein or free monomer) relative to the polymer chain, 
the particle is represented by an additional bead, not connected to the polymer chain, and 
subject only to hydrodynamic interactions between the particle and the chain. To avoid 
numerical instabilities in situations where the free particle overlaps with the polymer, in this 
case we use the Rotne-Prager-Yamakawa tensor (21), which modifies Eq. [3] for overlapping 
beads: 

y(r l3 ) = J(l- 9 Mr+ 3 ^^^} ifr y <2a. (6) 



v 32a/ 32a rf, 

i_ \ / lj _j 

Eq. [T]is discretized and integrated numerically using the Euler algorithm. The correlated 
stochastic contributions of Eq. H] are obtained from uncorrelated Gaussian noise by means 
of a Cholesky decomposition of the hydrodynamic matrix /x ^. In all the results below, 
lengths are measured in units of a, energies in units of k#T and time in units of r = 
a 2 /(kBTfio). The time step is At = 3 x 10~ 4 r, and a typical simulation lasts 10 9 steps, 
after an initial thermalization period of 10 6 — 10 7 steps. To reduce computational costs 
the Cholesky decomposition is only performed every 5 time steps. For a given chain length 
L = 2a(M — l) and persistence length l p , the quantities of interest are averaged over 25 — 2500 
different trajectories until the convergence is satisfactory. 



3 Theory 

3.1 Mean field theory for semiflexible polymer dynamics 

We briefly review the mean field theory for the dynamics of semiflexible polymers discussed 
in detail in Ref. (I22l ). In that earlier study the model was shown to provide an accurate 
description of internal polymer kinetics, validated through extensive comparisons with BD 
simulations. Moreover, the theory can be independently tested by comparison to recent 
fluorescence correlation spectroscopy experiments: without any fitting parameters, it exhibits 
excellent agreement for the mean square displacement (MSD) of tagged ends of single dsDNA 



fragments diffusing in solution (1231 l24j). 

The simplest description of a semiflexible polymer is the worm-like chain model: the 
polymer is represented by a continuous, different iable space curve r(s) of contour length L, 
where the arc- length variable s ranges from —L/2 to L/2. The associated elastic energy 
fAvLC, the continuum analogue of Eq. [5], is given by (1251 ): 
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The bending rigidity is k = Ipk^T ', and the tangent vector u = dr/ds is constrained by local 
inextensibility to unit length, u 2 (s) = 1 at each s. Since this constraint leads to nonlinear 
equations of motion, an alternative, approximate model is required. Within the mean field 



theory approach (1261 . 1271 ) the local constraint is relaxed and replaced by the global and end- 
point conditions dsu 2 (s)} = L and (w 2 (±L/2)) = 1 . The result is a Gaussian mean field 
Hamiltonian which incorporates a finite extensibility in addition to the bending term: 

CW[r(s)] = | f'^ds J^dsu 2 (s) + u (u 2 (L/2) + u 2 (-L/2)), (8) 

where e = 3l p k-BT/2, v = SksT/Alp, and v§ = 3ksT/A. In this form the Gaussian model 
exactly reproduces various lowest-order equilibrium averages of the worm-like chain, most 
importantly the tangent-tangent correlation function, and other derived quantities such as 
the mean square end-to-end distance. 

The dynamic theory for the Gaussian semiflexible polymer is based on the hydrodynamic 



pre-averaging approach of Ref. (I28l ) . analogous to that used for the Zimm model (I29l ) in the 
case of flexible chains. The time evolution of a point s on the polymer contour is governed 
by the Langevin equation: 



§- s r(s, t) = - /; /2 ds> y avg (s, s>) + i{ 8 , t), (g) 

(£0, t) ® t')) = 2k B T Vav g (s, s 1 ) 5{t - t'). 



Here we use the pre-averaged mobility tensor */Za Vg (s, s'), which is obtained from the standard 
Rotne-Prager tensor by averaging over all equilibrium configurations of the polymer. As seen 
in Eq. [3l the original Rotne-Prager mobility involves a dependence on the spatial distance 
between polymer points, and hence on the specific configuration of the chain. This would 
lead to nonlinear equations of motion, a problem which is resolved in the pre-averaging 
approximation, where the mobility depends only on the arc-length coordinates s and s'. For 
the mean field model of Eq. [H]the pre-averaged mobility has the form (j28l ): 



o r/ a e(|s-s'|-2a) / 6a 2 
2ajj, d{s — s ) H exp 



77 V /6tt 3 (t(|s - s'\) V 



(10) 



where a(l) = 2l p l — 2/ 2 (l — exp(— l/lpj)- The microscopic length scale a in the continuum 
theory corresponds to the monomer radius in the discrete BD simulations, and the unit step 
function in Eq. [TD] serves as a short- distance cutoff for the hydrodynamic interactions. 

The pre-averaged Langevin equation can be solved through a normal mode decomposi- 
tion, with the eigenmodes fulfilling free-end boundary conditions at s = ±L/2. This reduces 
Eq. [9] to a set of ordinary differential equations coupled by a hydrodynamic interaction ma- 
trix; the diagonalization of this matrix yields simple Langevin equations for the decoupled 
normal mode amplitudes P n (t) (with stochastic contributions Q n {t))\ 

-P (t) = Q (t), Ql P n{t) = -KPn{t) + Q n (t), U = 1,...,N, 

(Qm(t)Q mj (t')) = 2k B T8 ij 8{t-t')Q n 8 nm . 
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The vectors P n {t) and Q n (t) are related to r(s, t) and the stochastic velocities £(s, t) through 
the expansions r(s, t) = ^2 n=0 Pn(t)^n(s) and £(s, t) = J2n=oQn(t)^n(s), where the scalar 
functions \l/ n (s) are the decoupled normal modes. The modes are ordered in such a way that 



the eigenvalues A n (inverse relaxation times) increase with n. Following Ref. (1221 ). we set the 
high-frequency cutoff N for the mode number to N = \_L/8a\ + 1, which was shown to give 
good agreement at short times with BD simulations. At longer times, where the polymer 
fluctuations are at length scales much larger than the monomer radius a, the dynamics does 
not depend on the precise choice of the cutoff. A n and the fluctuation-dissipation parameters 
G ra can be directly derived from the tensor ^Jt^g evaluated numerically in the normal mode 
basis. Full details of this procedure, together with the explicit form of the normal modes 



ty n (s), are given in Ref. (1221 ) . 

Eq. [TT] for a given mode n can be mapped onto the well-known problem of a particle 
diffusing with friction constant ( in a harmonic potential of strength k centered at x = 



(I30l). We identify the mode amplitude P n (t) with the position x of the particle and set 
A n = k/(, n = This allows us to calculate the explicit form of the Green's function 
for the diffusive process, in other words the transition probability P(P n (t), P n (0); t) for 
observing the amplitude P n (t) after time t, starting from the initial amplitude P n (0): 

P(P n (t),P n (0);t) = (2nal(t)y 3/2 exp ^-^^^Tj . (12) 

The probability has the form of a spreading Gaussian, with time-dependent mean M n (t) 
and variance cr 2 (t) given by: 

M (t) = Po(0), M n (t) = P n (0)e- Ant , Vn > 0, (13) 
a 2 (t) = 2LD pol t, cr 2 n (t) = k B T6 n A- 1 (l-e- 2A » t ), \/n > 0. (14) 

Here D po \ = Q ^/ (s) 2 denotes the center-of-mass diffusion constant of the polymer coil. 
Modes with different mode number n evolve independently of each other. Transition prob- 
abilities between different polymer configurations, each of them corresponding to a certain 
(unique) set of normal mode amplitudes {P n (t)}, are therefore expressed as the product of 
the transition probabilities for the individual modes. As we are interested in the motion 
of a single point s on the polymer contour, the integrations over initial configurations with 
r(s, 0) = Vq and final configurations with r(s,t) = r can be readily performed. The result 
is a Gaussian Green's function specifying the conditional probability that point s reaches 
spatial position r in time t given that it starts at vq: 

G(r, r ; t) = (2nV(t)y 3/2 exp (- ^f ) , (15) 
V(t) =2D pol t + 2k B T 5^(1- e-^nis) 2 . (16) 



n=l n 



For t ^> A^ 1 the variance V(t) ~ 2D vo \t and hence the Brownian diffusion of the polymer's 
center-of-mass dominates the motion on large time scales. For smaller times, however, the 
contribution from internal polymer modes becomes important. The radial Green's function 
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G Ta d(r; t) for a particle starting at r = is obtained by integrating G(r, 0; t) over the surface 
of a sphere of radius r: 

Airr 2 ( r 2 \ 

G '- (r;t, = (^^ exp rw (17) 



In Sec. 14.11 we compare the time evolution of this transition probability for the case of the 
end-monomer (s = ±L/2) to BD simulation results. 

3.2 Hydrodynamic interactions between a target site and a free 
particle 

In the previous section we were able to describe the long-range hydrodynamic interactions 
between various points on a polymer coil, and their influence on the internal relaxation of 
the chain. A key simplifying feature in this analysis is the fact that the spatial relationship 
between any two points on the contour is constrained: their mean separation as the polymer 
fluctuates in equilibrium can be taken into account through the pre-averaging approximation. 
For the case of a free particle and a target site on the polymer coil, estimating hydrodynamic 
interactions is more difficult: the particle can drift away, and the strength of the interaction 
will be highly dependent on the initial conditions and the elapsed time. To understand the 
role of hydrodynamics in this situation, we first consider the simpler case of two freely diffus- 
ing spherical particles, before tackling the full problem of polymer-particle hydrodynamics. 
Though it may seem trivial, even the two particle case presents a challenging problem and 



can only be dealt with approximately (13ll . |32| ). 



We consider two non-interacting particles at positions 7*j(t), % = 1,2, described by the 
Langevin equations: 

(tM ® £,■(*')) = 2k B TS(t - 1') [%//iT + (i - <y V(r 12 ) n 

where m is the self-mobility of particle i, and hydrodynamics are expressed through the 
Rotne-Prager tensor \l (r 12 ), defined in Eq. [3X dependent on the inter-particle separation 
r 12 . Hydrodynamic interactions are long-ranged (oc r^ 1 ) and therefore negligible only at 
distances much larger than the sum of the particle radii, r 12 ^> ai + a 2 . For small separations, 
the stochastic motion of the particles is highly correlated, leading to a slowing down of their 
relative mobility. To get a realistic estimate of binding rates, where particles clearly have 
to approach each other, it is therefore necessary to take these hydrodynamic effects into 
consideration. 

The main quantity of interest is the radial Green's function for the relative motion of the 
particles, the probability that two particles starting from a distance ro = |ri 2 (0)| reach a 
distance r = |7"i 2 (t)| in time t. For comparison, we will consider the situation both with and 
without hydrodynamic interactions, labeling the respective Green's function G^ ad and G^ d . 
The second case, where the particles are totally decoupled and only the self-mobilities enter 
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into the stochastic correlations of Eq. CEBJ is trivial: 



GJU(r,r ;f) 



r oy /2nV*{t) 
V n (t)=2(D 1 + D 2 )t, 



2V n (t) 



g 2V n (t) 



(19) 
(20) 



where Di = fak^T is the diffusion constant of particle i. Note that the variance V n (t) is just 
1/3 of the MSD ((r 12 (t) - r 12 (0)) 2 ). 

In contrast, for the hydrodynamic case one cannot derive an analytical form for G^ ad . 
Thus we will have to resort to a heuristic approximation: we assume G^ ad has the same 
functional form as Eq. [19], but with a different variance V h (t), reflecting the slower relative 
motion of the particles. Since this variance is related to the MSD of ri 2 (t), we begin by 
evaluating this MSD. From Eq. [TH]and the definition of the Rotne-Prager tensor, Eq. [31 one 
can obtain the following exact relationship: 



((ria(t) - r 12 (0)) 2 ) = 6(D 1 + D 2 )t - 12 k B T 



d£ -—: 

o 6%r)r 12 {t') 



(21) 



=xt 



The second term on the right can be rewritten as — 12k B Txt, defining a time-dependent 
coupling parameter x which quantifies the slow-down in relative diffusion compared to the 
non-hydrodynamic case. This parameter involves both a time and ensemble average over 
the trajectory ri 2 (t), and hence is also dependent on the initial separation r . Our heuristic 
approach approximates x by an effective coupling function x(r ,t), where the time and 
ensemble averages have been replaced by a single ensemble average involving the Green's 



function G^ ad for the non-hydrodynamic system: 



X(r ,t) 



oo -I 

dr -Gl 
r 



n . ,n k vT 1 
- d(r ' ro;t)= 6^2^ 



erf | r ° ~ r ° 
6r \^/2VHt) / 



+ erf + 



(22) 



As in the pre-averaging approximation of the previous section, hydrodynamic interactions are 
cut off below the distance r a = a\ + a 2 where the particles overlap. This effective parameter 
shows the correct limiting behavior: x( r o ; ^) ^ t q 1 f° r short times where particles are close 
to their initial separation, and x( r o ; t) oc ((-Di + D 2 )t)~ 1 / 2 for long times where the particles 
have drifted far away from each other. In analogy to the non-hydrodynamic case, we define 
the variance V h (t) as 1/3 of the MSD, using the effective parameter x{ r o,t) instead of x- 
Thus the final form for our approximate hydrodynamic Green's function is: 



V h (t,r ) 



r ^/27iV h {t,r 
2(D l + D 2 -2 X (r ,t))t 



e 2Vll ( t .''o) — e 2V (*.n>) 



(23) 
(24) 



The role of the effective coupling parameter x is to reduce the relative mobility of the particles 
when they are near to each other and thus subject to strong hydrodynamic interactions; in 



DNA-Protein Binding Rates: Bending Fluctuations and Hydrodynamics 



9 



the long-time limit the hydro dynamic effects become negligible as the particles move to large 
separations. In Sec. 14.21 we test our approximation by comparing it to BD simulation results 
and show that it works well. 

In principle, the procedure outlined above to estimate \ can be iterated to produce 
higher-order approximations: in deriving x(r ,t) one can use the hydrodynamic Green's 
function G^ ad of Eq. [23] instead of G^ ad . This would lead to a better approximation for x, 
and hence a more accurate Green's function which could be input into the next level of 
the approximation. This iterative method converges quickly, so for simplicity we restrict 
ourselves to the first order results shown above. 

This approach for two freely diffusing particles can be generalized to the problem of a 
free particle and a polymer. The non-hydrodynamic case is again simple, with the relative 
motion of the particle and a point s on the chain described by a radial Green's function 
which has exactly the same form as Eq. [191 The only difference is the variance V n (t), which 
now includes the contribution of the polymer's internal modes (see Eq. IT6l) : 



Here D par = k&T /6irr)r pax is the diffusion constant of the free particle, and r par is the particle 
radius. 

For the hydrodynamic case, we can divide the complicated interactions between the 
particle and the chain into three parts: (i) the influence of the local region of radius a around 
the target site s on the free particle; (ii) the influence of the rest of the chain on the particle; 
(iii) the back influence of the particle on the entire polymer. For the specific problem we 
consider — association rates to a given target site — hydrodynamics plays a significant role only 
in the close vicinity of the target. Hence contribution (i) will dominate. The back influence 
in (iii) should be negligible for free particles comparable in size to the monomers in the 
chain, r par ~ 0(a), since the motion of the polymer is mainly governed by relaxation of the 
internal modes. The relative unimportance of (ii) is more subtle: one can take it into account 
in a more elaborate numerical evaluation of the pre-averaged MFT Langevin equations, 
but comparison to the simpler approximation discussed below does not show significant 
improvement with respect to BD simulations (which include all three contributions). Thus 
we can construct a simple estimate for the hydrodynamic Green's function by focusing on 
contribution (i). G^ ad has the same form as the two-particle case, Eq. [221 but with a variance 
given by: 



in analogy to Eq. [2H Here V n (t) is the non-hydrodynamic polymer-particle variance of 
Eq. [2SJ In this way we have accounted for both the internal fluctuations of the target 
site and, through x(ro, t), the decrease in relative mobility for the approaching particle. Our 
heuristic treatment of hydrodynamic effects is compared to BD simulation results in Sec. 14.31 




N-l 



(25) 



V h (t,r ) =V*(t)-4 X (r ,t)t 




(26) 
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3.3 Diffusion-controlled reaction rates 

In general, association rates are calculated by fin ding the steady-state solution of a diffusion 



equation with absorbing boundary conditions (1331. l34l). In our case, we can employ the simpler 



renewal approach (I35l ). which avoids the necessity of imposing absorbing boundaries: first 
passage times and binding rates can be directly extracted from the solution of the unbounded 
problem. This allows us to use the Green's functions derived in the last section for the 
relative motion of a target site on the polymer and a freely diffusing particle. The renewal 
approach is in principle equivalent to the Wilemski-Fixman formalism, and they share the 



same underlying approximations (1361 . 1371 ): (i) the Green's function for the time evolution 
of the particle-target distance is assumed to be that of a stationary Markov process; (ii) 
excluded volume between the reactands is ignored. We briefly review the main aspects of 
the approach below. 

For a stochastic process in one dimension, the Green's function g(x,x ;t) specifies the 
probability to find the particle at x at time t given a starting point x at time t — 0. If we are 
interested in paths that reach a boundary point x a in time t, the corresponding transition 
probability can be written convolution, 

g(x a ,x ;t)= dt' f(t';x a ,x )g(x a ,x a ;t-t'), (27) 
Jo 

where f(t;x a ,Xo) is the first passage time distribution: the probability of reaching x a in 
time t starting from xq at t = without passing through x a along the way. A Laplace 
transformation C acting on the time variable t leads to 

f{s;x a ,x ) = — <£> f{t;x a ,x ) = C [— — ], (28) 

gyx a ,x a ,sj g(x a ,x ai s) 

where Laplace transforms are denoted by £[/(£)] = f(s). Note that g(x, x$; t) is the Green's 
function in the absence of absorbing boundary conditions; we assume only that the proba- 
bility vanishes at x — > ±oo. 

For the various systems we consider, the binding process between two reactive particles 
is effectively one dimensional, depending only on the distance r between them. In the case 
of perfect absorption, particles bind when they collide for the first time — when the relative 
distance r reaches the absorption radius r a . As a result, the binding rate can be obtained 
from the first passage time distribution by integrating over all possible initial separations: 

POO 

k a {t) =4tt / dr r 2 f(t;r a ,r ). (29) 

J r a 

The steady-state rate k a reached at long times is extracted from the Laplace transform k a (s) 
using the final value theorem: 

k a = lim k a (t) = lim sk a (s). (30) 

For the simple case of two kinds of uncoupled Brownian particles with diffusion constants 
D\ and D 2 , the radial Green's function, Eqs. [1911201 yields: 

f(s;r a ,r ) = ^r ex P (-Oo ~ r ^\fj^)' ( 31 ) 
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with total diffusion constant D = D% + D<i- Applying Eq. [29j which corresponds to impose 
uniform initial concentrations at t — 0, one obtains the association rate: 

k a (t) = AnDr a (l + ~^=) , (32) 
V ynDt J 

which at long times reduces to the well-known Smoluchowski rate (jsl) ks = A.TiDr a con- 
stituting the upper limit for reaction rates governed by Brownian diffusion. For the more 
complicated polymer-particle case, the renewal method works analogously, but the Laplace 
transforms must be carried out numerically; the results are presented in section Sec. 14.41 



4 Results 

In order to validate the various analytical approaches described in Sec. El we test them against 
BD simulations. Since the binding rates depend sensitively on having good estimates for the 
Green's functions, we will focus on showing that the MFT approximation for the transition 
probabilities can reproduce the crucial physical effects: (i) the influence of internal polymer 
modes on the diffusion of the target site; (ii) the slow-down of relative motion between the 
free particle and target due to hydrodynamics. For simplicity, we concentrate in our analysis 
on one particular target site, the end-monomer of the chain. This has special relevance in 
biological processes, for example in the case of polymerization. 



4.1 Polymer motion 

We begin by considering just the internal relaxation of the polymer: the radial Green's 
function of Eq. [IT] for s = ±L/2 describing the diffusive motion of the polymer end-point. 
In Fig. [I] we compare this analytical MFT expression for G ra d(V, t) to histograms extracted 
from simulations of a chain with L = 100a, l p = 20a. The histograms are based on the 
analysis of 25 independent trajectories, each with 10 8 steps. The time evolution of the 
probability distribution is in excellent agreement with G Ta/ i(r,t) over time scales spanning 
four orders of magnitude. Note that the largest time-scale considered (t = 300 r) is less than 
the largest relaxation time of the polymer, T\ ~ 2 x 10 3 r, meaning that internal fluctuations 
dominate the polymer motion throu gho ut this entire time range. We illustrate the time 



evolution more directly in Movie SI (138l ). which depicts the motion of polymer end-points 
(red dots) superimposed from an ensemble of simulation trajectories, each starting at the 
origin (purple dot). The inset shows the evolving radial probability P(r,t) = G ra d(^; t)/r 2 
taken from the simulation histograms (red markers) compared to the MFT results (blue 
curve). Again we see a remarkable agreement, even in the tails of the distribution (which 
are visible on the logarithmic scale spanning seven decades of magnitude). Our comparison 



here is more detailed than in the earlier MFT study of Ref. (1221 ). since we consider the full 
transition probability and not just the end-point MSD. But the conclusion is the same: the 
MFT provides a highly accurate picture of internal polymer dynamics. 
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Figure 1: The probability G ra d( r , t) that the end-point of a polymer with length L = 100a 
and persistence length l p = 20a will diffuse a distance r in time t. The MFT predictions 
(lines) are compared to the results from BD simulations (black dots) for three different time 
regimes. Times are measured in units of r = a 2 / (k^T fi Q ) . 
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Figure 2: The probability G" ad (r, r ; £) that two spherical particles of radius a end up at 
separation r after time t given an initial separation r$ = 3a. For each time t, BD simulation 
data including hydrodynamics (black dots) are compared to two theoretical results: the non- 
hydrodynamic Green's function for decoupled particles labeled by a = n (Eqs. [H?l and |2"U|) . 
and the approximate hydrodynamic Green's function labeled by a = h (Eqs. I2"31 and . 
Times are measured in units of r = a 2 / {k^T hq) . 



4.2 Relative motion of two spherical particles including hydrody- 
namics 

Since the heuristic estimate for the hydrodynamic Green's function G^ ad (r, r ; t) of two freely 
diffusing particles (Eqs. |23H!M]) is the basis of our approach to polymer-particle interactions, 
we need to check that this estimate is reasonable. For this purpose, we performed BD 
simulations of two single spheres of radius a and extracted histograms from the variation of 
the inter-bead distance over time. In total 2400 independent trajectories were calculated, 
each starting at a distance r init = 3a and lasting 10 7 time steps. For given separations r 
and r, relevant transition events were identified, i.e. parts of the trajectories starting at a 
certain distance ro±oV/2 and ending at r±Sr/2, where 5r is the histogram binwidth. These 
were used to estimate the probability distribution corresponding to G^ ad (r, r ; t). In Figs. [2] 
and [3] we show the BD simulation results for two different r and various t together with 
G* ad (Eq. E3H2U) and the decoupled Green's function G* d (Eq. [191120]). The inhibition of 
relative motion due to hydrodynamics is clearly evident at short times: the distributions 
from the simulations are more narrowly peaked than the decoupled results, which spread 
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Figure 3: Same as in Fig. [2J but with initial separation r = 5a between the beads. 



out noticeably faster. The slowing down is correctly reproduced in the G^ ad curves, though 
it is slightly underestimated when compared to the simulation data. Despite the simplicity 
of the underlying approximation, G^ d still captures the essential features of the interaction. 
For longer times, when the inter-particle distance has become large and hydrodynamics plays 
a smaller role, G^ ad converges to G^ ad , and both agree with the BD results. 

4.3 Relative motion of a polymer and a spherical particle includ- 
ing hydrodynamics 

The final test for our Green's function approach concerns the approximate description of 
polymer-particle relative motion, contained in Eqs. [23]and[26j To extract the corresponding 
transition probabilities from BD simulations, we take a polymer with L = 100a, l p = 20a, 
allow it to thermalize, and then add a single free particle of radius a that is hydrodynamically 
coupled to the polymer beads. The particle is positioned at an initial distance r init = 3a 
from one of the polymer ends, and data is collected over 10 7 time steps. This procedure is 
repeated to obtain 2500 independent trajectories. As in Sec. I4.2[ relevant transition events 
are used to determine the probability distribution for diffusing from a separation ro to r 
in time t. The results, plotted in Figs. H] and [5], are qualitatively similar to the case of 
two particles: compared to the non-hydrodynamic case we see a reduced relative mobility 
between the target site and the particle for small separations. The hydrodynamic Green's 
function is quite close to the simulation results, again slightly underestimating the strength of 
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Figure 4: The probability G" ad (r, r ; £) that the end-monomer of a polymer with L = 100a, 
l p = 20a, and a spherical particle of radius a end up at separation r after time £ given an 
initial separation r = 3a. For each time t, BD simulation data including hydrodynamics 
(black dots) are compared to two theoretical results: the non-hydrodynamic Green's function 
labeled by a = n (Eqs. [T9l and 125]) . and the approximate hydrodynamic Green's function 
labeled by a = h (Eqs. [231 and 126]) . Times are measured in units of r = a 2 / '{k^T // ) ■ 
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Figure 5: Same as in Fig. HI but with initial separation r = 5a between the end-monomer 
and the free particle. 
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Figure 6: Schematic view of a bead-spring polymer with a reactive end and a free particle 
of radius r par . The particle is absorbed as soon as the separation r reaches the absorption 
radius r a = r par + a. 



the coupling. However, given the complexity of the correlations between the particle and the 
entire polymer coil, it is notable that we are able to get good quantitative agreement. Since 
the association rates, discussed in the next section, are derived directly from this Green's 
function, we should be able to obtain realistic estimates for the reaction process. 

4.4 Association rates 

The diffusion-limited reaction we consider is schematically illustrated in Fig. [6j free particles 
are absorbed by the reactive end-monomer of a chain as soon as the separation r reaches the 
absorption radius r a = r par + a. The association rate is calculated from the renewal approach 
in Sec. 13.31 with the only input being the radial Green's functions of Sec. 13.21 for the relative 
motion of the particle and the chain end. 

Since polymer fluctuations and the polymer-particle hydrodynamic coupling have com- 
peting effects on the association rate, it will be instructive to start with the simple case 
where the coupling has been turned off, i.e. we use the non-hydro dynamic Green's function 
G° ad of Eqs. [inland [23 The resulting association rate as a function of particle radius r par for 
a polymer with L = 1000a, l p = 50a is marked by red dots in Fig. [7J (For DNA, where a rs 1 
nm and l p ~ 50 nm, this would correspond to a strand of length 1 fim.) The Smoluchowski 
rate, 



which involves only the center-of-mass motion of the polymer coil and particle, is shown 
as a line; the corresponding numerical results of the renewal approach are shown as green 
triangles. The Smoluchowski value is the standard point of reference when considering DNA- 
protein interaction (l39l ). Clearly, internal polymer fluctuations have a significant impact, 
enhancing the association rates by 50 — 135% relative to ks over the range r par = 0.5a — 6a; 
the minimum at r par « 2a for the rates with internal polymer motion indicates that the 
effect of larger absorption radius r a = r par + a with increasing r par quickly dominates the one 
of decreasing diffusion constant _D par oc rZ^. 




(33) 
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Figure 7: Three estimates for the association rate of a free particle to the end-point of a 
polymer with L = 1000a, l p = 50a, as a function of particle radius r par . For DNA, where 
a ~ 1 nm and l p ~ 50 nm, the corresponding strand would have a length L = 1 /im. 
The rates are measured in units of k^T/6nn, or « 1.3 x 10 s M _1 s _1 for water at room 
temperature. Green triangles and solid line: Smoluchowski rate (Eq. I3"31 considering only 
the center-of-mass diffusion of the particle and the polymer coil, without hydrodynamic 
coupling. Red dots : the rates including the internal fluctuations of the polymer, without 
hydrodynamics (Eqs. [19] and [25]) . Blue crosses : the total rates including both fluctuations 
and hydrodynamics (Eqs. [T91 and I2"5]) . All the data marked by symbols is obtained from the 
numerical evaluation of the renewal approach (Eqs. |2"8][3"0"]) . 
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In the absence of hydrodynamics, we can use an earlier result of Berg (1l8l ) to do a 
consistency check on the renewal approach derivation of the association rates. For the case 
of a Gaussian Green's function, and assuming a Gaussian sink profile S(r) = exp (— 3r 2 /2r 2 ) 
instead of a perfectly absorbing boundary at r a , Berg derived an expression for the association 
rate in terms of the relative non-hydrodynamic variance V n (t): 



k Bcr g = (27r) 3/2 







dt (V»(t) + 2/3r 2 a ) 



-3/2 



(34) 



Plugging V n (t) from Eq. [25] into Eq. [3H we recover the non-hydrodynamic rates shown in 
Fig. [7] within a difference of 5% (comparison not shown). 

When hydrodynamic interactions are included together with internal polymer motion 
(Grad °f ^qs. E3] and [26] is used instead of G^ ad of Eqs. [[9] and l2"51) . the association rates 
are smaller than the non-hydrodynamic case, as seen in the data marked by blue crosses in 
Fig. [7J This is due to the inhibited mobility between the particle and the target at short 
distances. However the rate decrease is only ~ 15 — 20%, so the overall association rate is still 
30 — 100% larger than the Smoluchowski result for the range of particle sizes considered. The 
magnitude of the rate decrease is comparable to previous estimates derived for the simpler 
problem of a reaction between two spherical particles: Friedman obtained a 15% reduction 



due to hydrodynamics (13 ll ). while Deutch and Felderhof saw a 46% drop-off (j32l ). 

Since the conformational fluctuations depend on the mechanical properties of the chain, it 
is natural that the association rates will vary with the parameters L and l p that characterize 
the semiflexible polymer. In Fig. [8] we show contour diagrams of the association rate in 
terms of L/a and l p /a for two different particle radii: r par = a in the top panel, as would be 
the case in a polymerization process, and r par = 4a in the bottom panel, which corresponds 
to an average protein size for the case of DNA-protein interaction (with a = lnm). Fig. [9] 
depicts the same data, but in terms of percent rate increase over the Smoluchowski result 
ks of Eq. [33] Interestingly, while the absolute rate in Fig. [8] is almost independent of L for 
L ^> l p , it does show a slight maximum for l p « 20a when r par = a and l p ~ 30a when 
r par = 4a. This behavior as a function of l p is an intrinsic property of the polymer, because 
it appears also in the absence of hydrodynamic coupling (data not shown). The maximum 
is more prominent in the relative rates of Fig. [9] where the increase compared to ks rises 
sharply to ~ 30% (r par = a) and ~ 100% (r par = 4a) for long polymers with l p = 20a — 60a 
including double-stranded DNA {l p rs 50a). 



5 Discussion and summary 

In this paper, we have studied the association of free particles and a target site on a semi- 
flexible polymer — a class of diffusion-limited reactions of wide biological relevance, whether 
in polymerization of biomolecules or gene regulation by DNA-binding proteins. We focused 
on two competing effects that are, with the lone exception of the segmental diffusion con- 
sidered in Ref. (1181 ). entirely neglected in existing theories for these processes: the bending 
fluctuations of the polymer in equilibrium, which enhance the association rate, and the 
hydrodynamics between the polymer and particle, which reduces it. Quantifying the fluc- 
tuations required an accurate description of internal polymer motion, available through the 
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Figure 8: The association rate k a (Eqs. [28ll30l with Eqs. [23] and [26]) between a free particle 
and a polymer end as it varies with the polymer contour length L and persistence length 
l p . Top panel: for particle radius r par = a. Bottom panel: for r par = 4a. The rates include 
the effects of polymer fluctuations and hydrodynamic coupling, and are measured in units 
of k-^T /QTcn, or w 1.3 X 10 8 M _1 s _1 in water at room temperature. 
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Figure 9: Same as in Fig. [HJ but showing the percent rate increase over the Smoluchowski 
result, ks (Eq. [33]), instead of the absolute rate. 
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Figure 10: Schematic representation of two scenarios for a bimolecular reaction: I) both 
species are free in solution; II) one of the species is attached to the end of a semiflexible 
polymer. The comparison of the diffusion-limited association rates k\ and fc* 1 should highlight 
the role of internal polymer motion in such reactions. 



mean field theory of semiflexible polymer dynamics. For the hydrodynamics, we developed 
a simple heuristic estimate to model the decrease in mobility when two diffusing objects ap- 
proach each other. The end result of the competition sensitively depends on the mechanical 
properties of the polymer and the size of the reactive particle: we see a maximal increase 
over the Smoluchowski rate for persistence lengths l p = 20a — 60a; this increase is ~ 30% for 
the case of small particles, and ~ 100% for particle sizes typical of regulatory proteins. As a 
function of particle size, the reaction rate displays a minimum at a particle radius equal to 
the polymer diameter; for larger proteins the probability of hitting the DNA increases. Since 
most measured DNA-protein reaction rates are on the order of the Smoluchowski rate ([]]), 
the effects described here will be crucial to getting a complete quantitative theory of these 
reactions. We have neglected the possibility of encounters between the protein and the target 
site on the DNA that do not immediately lead to a reaction, for example due to orientational 
constraints. This is clearly an oversimplification, which however can be easily rectified by a 
reaction factor smaller than unity (40). 

Directly testing our theoretical predictions may be possible with existing experimental 
techniques. In particular, one can study the role of polymer fluctuations on association 
rates by comparing two different setups involving bimolecular diffusion-limited reactions, 
shown schematically in Fig. [TUJ In scenario I both reactands are free in solution, while 
in scenario II one of the reactands is bound to the end of a semiflexible polymer. The 
corresponding rates k\ and fc„ could be measured through simple kinetic experiments or 
using fluorescence microscopy. In the latter case, the reactands would be fluorophores and 
quencher molecules, and the rate of the diffusion-limited quenching reaction can be extracted 
from the decay of the total fluorescence signal over time. As shown in Table [TJ our theory 
predicts a ratio k l }jk\ ~ 0.77 for the two scenarios using r par = a and a 1 urn strand of 
DNA. This is clearly distinguishable from the Smoluchowski prediction, k l ^/k l a = 0.5, which 
ignores hydrodynamics and internal fluctuations. 
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Rate estimates including 


k 1 

K a 


*? 




center of mass motion 


50.27 


25.14 


0.5 


plus internal polymer motion 




39.83 




plus hydrodynamics 


42.52 


32.71 


0.77 



Table 1: Estimates for diffusion-controlled bimolecular association rates in the experimental 
scenarios shown in Fig. [TTJ The polymer is assumed to be a 1 rim strand of dsDNA, and 
the particle radius r par = 1 nm. The rates k a are given in units of h^T /6Trr], or ~ 1.3 x 10 s 
M~ 1 s _1 in water at room temperature. 
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